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Abstract 


We  present  a  new  class  of  algorithms  for  determining  whether 
there  exists  a  point  x  e  Rn  satisfying  the  nonlinear  inequality  con¬ 
straints  c. (x)  <0,  i  =  1,  . ..,  m,  subject  perhaps  to  satisfying 
linear  inequality  constraints  £j(x)  <  0,  j  =  l,...,k  which  are  known 
to  be  feasible.  Our  algorithm  consists  of  solving  a  sequence  of 
linearly  constrained  optimization  problems,  using  a  sequence  of 
objective  functions  4>(x,p)  which  are  at  least  twice  continuously 
differentiable,  and  which  are  generated  by  monotonically  increasing 
the  value  of  the  non-negative  parameter  p.  It  is  shown  that  in 
almost  all  cases,  once  p  reaches  or  exceeds  some  finite  value,  that 
the  solution  to  the  linearly  constrained  optimization  problem 
either  is  a  feasible  point,  or  establishes  the  infeasibility  of  the 
set  of  constraints.  Computational  results  are  presented  in  which 
the  algorithm  performs  satisfactorily  on  feasible  and  on  infeasible 
systems . 
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1.  Introduction 

In  this  paper,  we  present  and  analyze  a  new  class  of  algorithms 
for  determining  wither  a  system  of  linear  and  nonlinear  inequality 
constraints  has  a  feasible  point.  The  general  problem  we  consider  is 

Given  D  e  Rn,  c- :  D  +  R,  i  *1 .  in 

Does  there  exist  any  x  e  such  that  (1.1) 

c.j (x)  <0,  i  =1,  ....  m? 

The  m  functions  c^{x)  are  called  constraint  functions,  or  simply  con¬ 
straints.  If  the  answer  to  (1.1)  is  yes,  the  set  of  constraints  {c^} 
is  said  to  be  feasible,  otherwise  it  is  said  to  be  infeasible.  These 
definitions  are  made  formally  in  Definition  1.1. 

Problem  (1.1)  has  many  applications.  Our  original  interest 
came  from  the  static  software  validation  project  at  the  University  of 
Colorado  at  Boulder,  where  (1.1)  arises  when  determining  satisfiability 
of  path  constraints  (see  e.g.,  [3]).  In  this  case,  many  of  the  con¬ 
straints  are  linear.  Since  problem  (1.1)  is  easily  solved  in  the 
linear  case,  it  is  advantageous  to  treat  linear  constraints  separately 
from  nonlinear  ones.  Also,  in  most  applications  of  (1.1)  it  is 
necessary  that  the  algorithm  actually  return  a  feasible  x  when  the 
answer  to  (1.1)  is  yes.  Thus  the  form  of  problem  (1.1)  we  will  address 
is: 

Given  D  e  Rn, 

nonlinear  constraints  c..  :  D  -*■  R,  i  =  1,  . . . ,  m 
linear  constraints  l- :  D-*R,  j  =  l,  ...,k 

J 

I)  Determine  whether  there  exist  any  x  e  D  such  that 
tj(x)  s  0,  j  =  1,  ..., 

(and  if  so:) 


k 


(1.2) 
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m 

,  k 


Part  I  of  problem  (1.2)  is  simply  phase  I  of  linear  programming, 
and  there  exist  many  reliable  codes  to  solve  it  (see  e.g.,  [4]).  Thus 
o ur  research  is  concerned  solely  with  the  solution  of  part  II  of  (1.2) 
in  the  case  when  the  linear  constraints  are  feasible.  It  should  be 
noted  that  our  theoretical  results  remain  valid  and  our  computer  al¬ 
gorithm  performs  just  as  well  in  the  case  k  =  0,  i.e.,  no  linear  con¬ 
straints  . 

The  reader  may  notice  from  the  wording  of  problem  (1.2)  that  we 
have  attempted  to  construct  a  global  algorithm  i.e.,  one  which  solves 
(1.2)  over  the  entire  domain  D,  as  opposed  to  a  local  algorithm,  i.e, 
one  which  solves  (1.2)  in  a  neighborhood  of  the  starting  point.  While 
we  realize  v.nas  a  foolproof  global  algorithm  for  sucn  a  nonlinear  pro- 
d-rr.  is  i  :p:-  s;  we  have  developed  our  algorithm  with  a  global 
solution  n  mind.,  and  in  practice  have  been  rather  successful.  We  have 
alv  assumed  v  rttia l  !_v  no  special  structure,  inducing  concavity  or 
convexity,  c it  the  nonlinear  constraints  c^(x).  In  the  applications  we 
have  seen,  it  is  reasonable  to  assume  that  each  c^  is  twice  continuously 
d i  f f  m-i-si.  i rd . * Ha  s,  che  computer  program  described  in  this  paper 
makes  use  of  t:  i ■  assumption.  However,  the  theoretical  results  of  this 
paper  do  not  even  rrginre  c^  continuous.  Other  than  this  we  make  no 
assumption  about  the  nonlinear  constraints. 

There  does  not  appear  to  be  much  previous  work  on  problem  (1.2), 
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especially  in  the  generality  in  which  we  are  interested.  A  paper  of 
Robinson  [  n]  discusses  a  local  method  for  finding  a  feasible  point, 
assuming  one  starts  close  enough  to  one.  Several  other  people,  includ¬ 
ing  Bertsekas  [1]  and  Charalambous  and  Conn  [2]  have  developed  algorithms 
wilier,  can  be  applied  to  problem  (1.2).  The  first  approach  would  require 
optimization  subject  to  nonlinear  inequality  constraints,  the  second, 
optimization  of  an  objective  function  which  is  not  continuously  differen¬ 
tiable.  Our  approach  differs  strongly  in  that  its  optimization  problems 
involve  twice  continuously  differentiable  objective  functions,  and  only 
linear  constraints. 

The  remainder  of  this  paper  is  organized  as  follows.  In  Section  2 
■o  m.rnvnue  a  no  present  a  new  class  of  algorithms  for  solving  problem 
! L . 2 )  cart  11.  It  involves  a  "weighting  function"  which  is  only  required 
to  -i L  i sfy  some  very  general  properties.  In  Section  3  we  prove  a  num¬ 
ber  :>f  theoretical  properties  of  our  class  of  algorithms,  which  show 
that  they  work  successfully  in  almost  all  cases.  Section  4  discusses 
i  e  c. «> ice  of  the  weighting  function  in  practice.  Section  5  describes 
c 'specis  of  our  computer  algorithm,  and  presents  some  computational 
res  Jr  ..  In  Section  6  we  summarize  our  work,  and  discuss  the  directions 
<’i  one  continuing  research  in  this  area. 

'Jo  conclude  this  section  with  a  formal  definition  of  the  con- 
opus  of  feasibility  and  infeasibility,  which  were  described  in  the 
ope’i i in;  paragraph.  In  addition,  we  introduce  the  terms  strictly  feas¬ 
ible.  strictly  infeasible  and  exactly  feasible,  which  will  be  important 
n  our  discussions  in  Sections  3  and  4.  Informally,  a  set  of  constraints 
is  strictly  feasible  if  it  has  a  feasible  point  where  in  addition  no 
constraint  is  equal  to  zero;  it  is  strictly  infeasible  if  it  is 
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infeasable  and  dots  not  come  arbitrarily  close  to  feasibility.  It  is 
exactly  feasible  if  it  is  feasible  but  not  strictly  feasible. 

Definition  l.i.  wot  D  e  Kn,  anc  let  c-  :  D  -+  R,  i  =  1,  . ..,  m.  The 
Set  of  const;  a  i fits  ■  is  said  to  be 

feasible  V  tr.er..-  exists  x  c  D  such  that  c.(x)  <  0,  i  =  1,  . ..,  m 
infeasible  if  -c.-  is  not  feasible,  i.e.,  for  all  x  e  D, 

max  { c .  ;  ' •  0 

i 

1  i  •:  n 

strictly  f>.-  i  •  i > «t  if  there  exists  e  <  0,  x  t  D  such  that  c^x)  <  e, 

strictly  inf--'.  mt  ie  if  there  exists  6>0  such  that  for  all 

x  i  0 ,  ..tax  •  a  .• 1  /  ■ 
i  ;  i 

exact"  v  i.asjb’p.  if  .  c.i  is  feus i ole  but  not  strictly  feasible. 
Given  an  .  0 ,  o-  is  said  tc  be 

strictly  .  asioje  if  tnore  exists  x  c  D  such  that  c.(x)  s  e, 

•  ‘  1  ,  •  •  •  ■  !■ 

(i i  '/■.  ii  a  ■  u  'tiiu  u'  .j a 

wi  y  m  ■  f  feasibility  ii  tnere  exists  x  c  D  such  that 

c  . '  a  )  .  i  i  .  ....  'ii. 

1 

dr. '.a  that  if  tne  domain  is  closed  and  bounded,  then  a  set  of 
constraints  wrier-  is  infeasiblf  must  be  strictly  infeasible.  On  a 
computer,  it  i  reasonable  to  assume  that  the  domain  is  closed,  with 
ra  rosnone'  c  o<  x  bounded  in  abso.  ite  value  by  the  largest  floating 
flirt.  number  e  -  *iie  rnarnine.  : nus ,  it  is  reasonable  to  assume  in 
ict  ice  ■  ha t  <; '  i  i  i feasible  problems  are  strictly  infeasible. 
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In  this  section,  we  introduce  our  new  class  of  algorithms  for 
solving  problem  (1.2),  part  II  (henceforth  called  (1.2II)).  Because 
it  is  somewhat  complex,  we  motivate  it  by  first  describing  an  algorithm 
for  a  simpler  problem.  Our  final  algorithm  is  built  around  this  simpler 
one. 

Suppose  first  we  wanted  to  solve  problem  (1.21 I)  with  only  one 
nonlinear  constraint,  i.e.,  m  =  1.  One  can  do  this  using  well-known 
techniques,  by  simply  solving  the  linearly  constrained  optimization 
problem 

min  c^x) 

xeD  (2.i) 

subject  to  £.(x)  <0,  j  =  1,  ....  k 

J 

If  at  the  minimum  x*  of  (2.1),  Cj(x*)  is  nonpositive,  then  the  set  of 
constraints  {c^,^,...,^}  is  feasible,  otherwise  it  is  infeasible. 
Notice,  however,  that  if  we  desire  a  global  answer  to  the  feasibility 
problem  (i.e.,  over  all  xeD),  then  we  will  require  the  global  minimum 
to  (2.1).  Also,  if  c^(x)  is  not  convex,  then  our  linearly  constrained 
optimization  routine  cannot  assume  convexity. 

A  computer  program  to  solve  (2.1)  under  these  conditions  is, 
perhaps  surprisingly,  the  key  piece  of  software  in  our  algorithm  to 
solve  the  more  general  feasibility  problem  ( 1 . 2 1 1 ) .  This  is  because 
our  algorithm  for  solving  ( 1 .21 I )  turns  out  to  consist  of  solving  a 
series  of  problems  of  form  (2.1).  Furthermore,  if  we  can  successfully 
find  the  global  minimum  to  the  problems  of  form  (2.1),  then  it  turns 
out  that  we  can  find  the  global  solution  to  (1.2II).  (Otherwise  we 
at  least  solve  ( 1 .21 1 )  over  whatever  domain  our  linearly  constrained 
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optimization  algorithm  is  actually  able  to  find  its  minima  on.)  Thus, 
the  global  aspect  of  our  original  problem  is  reduced  to  a  somewhat 
more  manageable  framework. 

Our  algorithm  for  solving  (2.1)  is  discussed  in  Section  5.  In 
brief,  we  use  the  best  existing  techniques  for  finding  a  local  minimum 
to  (2.1)  (see  e.g.,  [6,7,8]  )  ,  modified  if  necessary  in  the  face  of 
nonconvexity,  and  then  restart  in  some  other  region  of  the  domain  if 
either  curvature  information  accumulated  during  the  local  algorithm,  or 
a  sampling  of  points,  indicate  that  a  lower  minimum  may  be  found  there. 
The  routine  has  been  quite  successful  in  our  tests. 

Now  consider  the  general  case  of  problem  (1.2II)  with  m  >  1.  It 
can  also  be  reduced  to  the  solution  of  a  linearly  constrained  optimi¬ 
zation  problem  in  at  least  two  ways.  While  they  are  not  quite  what  we 


want,  they  lead  the  way  to  our  algorithm.  The  first  way  is  to  solve 


min 
x  c  D 


m  . 

♦l(x)  =  X  cj(*) 


(2.2) 


subject  to  £-(x)  <  0,  j  = 1,  . . . ,  k 

J 


subject  to  i.(x)  s  0,  j  *  1,  ...,  k. 

J 

If  the  minimum  value  of  <p2  is  nonpositive,  then  (1.2II)  is  feasible, 
otherwise  it  is  infeasible. 
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The  problem  with  both  of  these  approaches  is  that  they  involve 
objective  functions  <j>j  and  which  are  not  continuously  differentiable. 
Standard  software  for  linearly  constrained  optimization,  including  ours, 
does  not  work  on  such  problems,  and  modifying  it  to  be  successful  is 
difficult.  However,  (2.2)  especially  suggests  a  related  idea  which  is 

the  basis  of  our  algorithm:  to  solve  a  problem  of  the  form 

m 

min  *(x)  =  l  w(c.(x)) 

xeD  1  =  1  (2.4) 

subject  to  lAx)  <  0,  i  =  1,  ....  k 

J 

where  w:  R  R  is  a  weighting  function  which  penalizes  positive  con¬ 
straint  values  and  rewards  negative  constraint  values.  The  type  of 
weighting  function  we  will  use  is  drawn  in  Figure  2.1,  and  formally 
defined  in  Definition  2.1.  The  most  obvious  example  is  w(y)  =  ey-l, 
but  many  other  possibilities  are  also  discussed  in  Section  4. 


Figure  2.1  -  A  weighting  function 

Definition  2.1.  A  weighting  function  is  a  function  w  :  R  -*•  R  which  obeys 

1)  w  is  twice  continuously  differentiable. 

2)  w(0)  »  0. 

3)  w  is  monotonically  increasing  (i.e.,  w'(y)  >  0  for  all  ye  R). 


4)  w  is  strictly  convex  (i.e.,  w"  (y)  >  0  for  all  y  e  R). 

5)  there  exists  finite  l  <  0  such  that  lim  w(y)  =  l. 

y  “ 

6)  lim  w(y)  =  +  ® 

y  -►  +  oo 

Note  that  property  6  is  actually  implied  by  the  other  properties. 

Since  w  is  twice  continuously  differentiable,  problem  (2.4)  can 
be  solved  using  our  linearly  constrained  optimization  software.  Now 
let  us  see  what  information  its  solution  gives  us  about  the  original 
feasibility  problem  (1.2II).  Let  us  call  the  point  at  which  (2.4) 
attains  its  minimum  x*.  If  <j>(x*)  >  0  then,  just  as  in  problems  (2.2) 
and  (2.3),  (1.2II)  must  be  infeasible.  This  is  simply  because  if  a 
feasible  point  x  exists,  then  (j>(x)  <  0  by  properties  2  and  3  of  w,  and 
so  <f>(x*)  would  have  to  be  negative.  However,  if  4>(x*)  <  0  then  there 
are  two  possibilities:  x*  may  be  feasible,  or  some  constraints  may 
be  positive  and  others  negative  at  x*.  In  the  final  case,  the  solu¬ 
tion  to  (2.4)  does  not  give  us  the  solution  to  feasibility  problem 
(1.2II). 

These  possibilities  are  summarized  below: 

Let  x*  be  the  solution  to  (2.4)  where  w:  R  -+  R  obeys  Definition 
2.1.  Then  exactly  one  of  the  following  is  true: 

1.  4>(x*)  >  0  :  this  implies  (1.2II)  is  infeasible. 

2.  ${x*)  s  0  and  x*  is  feasible:  a  feasible  point  to  (1.2II) 
has  been  found. 

3.  4>(xJ  s  0  and  x*  is  infeasible:  (1.2II)  may  be  feasible  or 
infeasible. 

Only  the  third  case  is  unsatisfactory.  To  remedy  it,  we  make 
one  final  addition  to  our  process:  we  solve  a  sequence  of  problems 
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of  form  (2.4),  using  a  sequence  of  weighting  functions  obeying 
Definition  2.1.  These  weighting  functions  successively  penalize 
constraint  violation  more  and  more,  and  reward  constraint  satisfaction 
less  and  less.  In  this  manner,  they  attempt  to  make  any  feasible  point 
more  attractive  to  problem  (2.4)  than  any  infeasible  point.  A  rep¬ 
resentative  sequence  of  such  weighting  functions  is  drawn  in  Figure  2.2; 
it  is  shown  how  to  easily  construct  such  sequences  in  Lemma  2.2.  The 
key  results  of  this  paper  (in  Section  3)  then  show  that  by  solving 
such  a  sequence  of  problems  (2.4),  we  eventually  end  up  with  a  solu¬ 
tion  x*  to  (2.4)  which  falls  into  Case  1  or  2  above,  instead  of  Case  3. 
Thus,  we  solve  the  feasibility  problem  ( 1 . 21 1 ) .  (In  practice,  at 
least  over  the  domain  on  which  we  solve  the  minimization  problems.) 
Intuitively,  this  may  also  be  viewed  as  saying  that  as  our  weighting 
function  becomes  shaped  more  like  a  backward  "L,"  problem  (2.4) 
comes  to  resemble  problems  (2.2)  and  (2.3). 


A  family  of  weighting  functions  like  the  one  in  Figure  2.2  is 
constructed  using  a  standard  result  from  convexity  theory  (  see  e.g., 
[12]):  if  w(y)  is  any  function  satisfying  Definition  2.1,  and  { p^ } 
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is  any  monotom'cally  increasing  sequence  of  positive  real  numbers, 
then  the  functions 

w.j(y)  =  w(pi  .  y)  i  =  1,  2,  ...  (2.5) 

have  monotonical ly  increasing  values  on  each  nonzero  y  (as  well  as  each 
trivially  satisfying  Definition  2.1).  This  is  proven  in  somewhat 
more  generality  in  Lemma  2.2. 


Lamia  2.2  Let  w:  R  -*•  R  be  continuously  differentiable  and  strictly 
convex,  and  assume  w(0)  =  0.  Let  p+  >  P  >  0.  Then  for  any  y  e  R, 
w(p+-y)  w(p*y) 

"  P+  "  " 

with  equality  if  and  only  if  y  =  0. 

Proof:  If  y  =  0  the  lemma  is  trivially  true.  Now  suppose  y  *  0. 

Since  w  is  continuously  differentiable  and  strictly  convex, 

w(z)  -  w(x)  >  w‘(x)  (z-x)  (2.6) 

for  any  *,zeR,  x*z.  Substituting  x  =  p-y,  z  =  p+*y  into  (2.6)  gives 
w(p+y)  -  w(py)  >  w'(py)  •  (p+-p)y.  (2.7) 

Also,  substituting  x  =  py,  z  =  0  into  (2.6)  and  using  w(0)  =  0  gives 
-w(py)  >  -w' (py)  •  py 


or 


w'(py)-y  2  w(py) 

P 

Combining  (2.7)  and  (2.8) 

w(p+y)  -  w(py)  >  w(py)  (p+-p) 
P 

or  equivalently 


(2.8) 


.  “(WIP*  -  w(py) 


W^P+-'1^  >  w(py)  ■ 
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By  solving  a  sequence  of  problems  of  form  (2.4),  with  the  weight¬ 
ing  functions  constructed  according  to  (2.5),  we  have  completed  the 
development  of  our  class  of  algorithms  for  solving  feasibility  pro¬ 
blem  ( 1 .21 I ) .  It  is  given  below.  The  parameter  p  is  called  by 
that  letter  because  it  resembles  the  penalty  parameter  found  in 
various  mathematical  programming  algorithms.  D,  c-(x)  and  £.(x) 

*  J 

are  defined  as  in  problem  (1.2),  and  we  assume  that  the  linear  con¬ 
straints  have  been  found  to  be  feasible  in  problem  (1.2),  part  I. 

Algorithm  2.3  -  -  A  weighting  function-penalty  algorithm  for  de¬ 
termining  feasibility  of  a  set  of  nonlinear  and 
linear  inequality  constraints  (problem  ( 1 . 21 1 ) ) . 

1.  p  :  =  0 

1  m 

2.  solve  min  <j>(x,p)  =  -  l  w(p-c,(x)) 

r  i  =  1  1 

x  e  D  L 

subject  to  4-(x)  s  0  j  =  1,  ...»  k 

J 

say  the  minimum  occurs  at  x*(p) 

3.  if  4>(x*(p),p)  >  0,  stop  (constraints  are  infeasible) 
else  if  x*(p)  feasible,  stop  (feasible  point  found) 

else  increase  p  and  return  to  step  2 


Note  that  Algorithm  2.3  starts  with  p  =  0.  Since  w(0)  =  0, 
4> (x ,0)  =  0/0,  but  L'Hopital's  rule  shows  that 


<t>(x»0)  =  w'(0)  •  l  c^x) 
i=l 


(2.9) 


where  w'(0)  is  some  nonzero  constant.  So  (2.9)  is  the  first  objective 
function  used  in  step  2  of  Algorithm  2.3.  There  are  two  advantages 
to  starting  with  this  objective  function.  First,  the  first  linearly 
constrained  optimization  problem  to  be  solved  becomes  a  reasonably 
nice  one,  and  yet  its  solution  often  establishes  feasibility  or 
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infeasibility  in  practice.  Second,  when  m  »  1,  Algorithm  2.3  reduces 
to  the  algorithm  proposed  at  the  beginning  of  this  section  for  solving 
(1.2II)  with  one  nonlinear  constraint. 

The  most  important  fact  about  Algorithm  2.3  is  that  it  will  now 
be  shown  in  Section  3  that  it  is  guaranteed  to  produce  a  feasible  x, 
or  establish  infeasibility  of  the  nonlinear  constraints  as  soon  as  p 
reaches  or  exceeds  some  finite  threshhold.  (There  is  one  exceptional 
case.)  Thus,  while  one  might  fear  that  Algorithm  2.3  would  only  give 
the  answer  to  the  original  feasibility  problem  ( 1 . 21 1 )  in  the  limit 
as  p  -+  <*>,  it  in  fact  gives  the  answer  for  finite  p,  and  this  is  the 
key  reason  why  it  is  an  effective  computational  algorithm.  We  there¬ 
fore  feel  that  our  class  of  algorithms  is  much  more  similar  in  its 
penalty  characteristics  to  the  "augmented  Lagrangian"  or  "method  of 
multipliers"  methods  for  constrained  optimization  (see  e.g.,  [10]), 
which  are  also  guaranteed  to  terminate  satisfactorily  for  finite  p, 
than  to  the  classical  penalty  or  barrier  algorithms  for  constrained 
optimization  (see  e.g.,  [  5,9])  which  only  produce  the  correct 
answer  in  the  limit  as  p  -*•  ».  It  will  be  seen  in  Section  5  that 
typically  few  (1  to  3)  values  of  p  are  required  in  practice. 
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3.  Termination  properties  of  the  new  class  of  algorithms 

In  this  section,  we  establish  the  key  property  of  the  class  of 
algorithms  introduced  in  Section  2  (Algorithm  2.3)  for  determining 
feasibility  of  a  set  of  nonlinear  inequality  constraints:  that  in 
almost  all  cases,  the  algorithm  is  guaranteed  to  either  produce  a 
feasible  point,  or  establish  infeasibility,  as  soon  as  the  parameter 
p  reaches  or  exceeds  some  finite  value.  This  is  proven  in  Theorems 
3.4  -  3.6.  Before  this  we  prove  two  fundamental  properties  of  our 
class  of  algorithms.  We  conclude  the  section  by  discussing  the  termin¬ 
ation  criteria  which  are  used  in  the  computer  implementation  of  our 
algorithm. 

For  the  purposes  of  this  and  the  next  section,  there  is  no  reason 
to  make  explicit  mention  of  the  linear  constraints  £.(x),  j  =  1,  ....  k. 

J 

Rather,  we  will  assume  that  the  domain  D  has  been  restricted,  if 

necessary,  to  exclude  any  point  x  where  some  £,(x)  >  0.  This  simplifies 

J 

our  notation,  and  allows  us  to  use  the  definitions  of  feasibility  and 
infeasibility  from  Definition  1.1  without  making  any  modification  due 
to  the  linear  constraints.  Using  this  convention,  an  iteration  of 
Algorithm  2.3  consists  simply  of  minimizing  $>(x,p)  over  xeD.  We 
will  again  define  the  point  when  this  minimum  occurs  as  x*(p),  and 
for  brevity  will  refer  to  the  minimum  value  <f»(x*(p),p)  as  <j,*(p). 

This  notation  is  summarized  below. 

Definition  3.1  Let  D  e  Rn,  c.  :  D  -*•  R,  i  *  1 ,  . . . ,  m,  and  let 
w  :  R  R.  For  p  >  0,  the  composite  weighting  function  $(x,p) 
is  defined  by 

i  m 

♦(x.p)  *  -  l  w(p-ci(x)). 
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For  p  s  0,  $(x,p)  is  defined  by 

m 

4>(x ,0)  =  lim  4(x,p)  =  w'(0)  l  c.(x). 
p-*-0  i=l  1 

For  p  ^  0,  we  define 

x*(p)  =  min  <|){x,p) 

X  e  D 

and 

t>*(P)  =  <fr(x*(p),p). 

In  Lemma  3.2  we  prove  a  property  of  our  class  of  algorithms 
which  was  informally  indicated  in  Section  2:  that  if  (c^ }  is  feasible, 
then  <t>*(p)  is  nonpositive  for  all  p  2  0.  This  leads  directly  to  the 
criterion  which  is  used  to  establish  infeasibility  in  our  class  of 
algorithms:  if  4>*(p)  is  positive  for  any  p  >  0,  then  { c^ }  must  be 
infeasible. 

Note  that  the  lemmas  and  theorems  of  this  section  assume  only 
those  properties  of  w(y)  which  they  actually  use.  However,  taken 
together  they  require  all  the  properties  from  Definition  2.1  except 
for  w"  (y)  continuous.  This  is  required  only  by  the  linearly  con¬ 
strained  optimization  algorithm.  Recall  also  that  while  in  theory 
our  results  are  true  over  the  entire  domain  xeD,  in  practice  they 
may  only  be  valid  over  whatever  domain  D  5  D  the  linearly  constrained 
optimization  routine  can  in  fact  find  the  minima  of  4 ( x ,p ) .  That  is, 
if  the  algorithm  finds  a  feasible  point  the  problem  is  of  course 
feasible,  but  a  claim  of  infeasibility  may  only  be  valid  over  the 
domain  on  which  we  are  actually  able  to  solve  the  minimization  problems. 

Lenina  3.2  Let  w  :  R  -*■  R  satisfy  w(0)  =  0  and  w(y)  <  0  for  all 
y  <  0,  and  let  c^(x),  <j>(x,p),  x*(p),  <>*(p)  be  defined  by  Definition 
3.1.  If  {c.}  is  feasible,  then  for  any  p  >  0,  <t>*(p)  <  0.  If  {c..}  is 
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strictly  feasible,  then  for  any  p  0,  $*(p)  <  0. 

Proof:  If  (c.)  is  feasible,  then  there  exists  some  x  such  that 
cn.(x)  s  0,  i  =  1,  m.  Thus,  4>(x,p)  <  0  for  any  p  >  0,  and  since 
by  the  definition  of  <j>*(p),  <t>*(p)  <  4>(x,p),  we  have  cj.*(p)  <  0.  The 
strictly  feasible  case  is  proven  identically. 

A  question  which  might  arise  following  Lemma  3.2  is:  if 
<f>*(p)  >  0  for  sane  fixed  p  >  0,  is  it  possible  that  for  some  p+  >  p, 
4>*(p+)  ^  0?  If  this  were  possible,  then  our  test  for  infeasibility 
would  be  rather  chancy,  as  we  might  "skip  by"  the  choice  p  which  would 
establish  infeasibility,  and  use  instead  the  choice  p+  which  doesn't. 
Luckily  this  is  not  possible,  because  it  is  an  easy  consequence  of  the 
relation  of  successive  weighting  functions  (Lemma  2.2)  that  <j>*(p) 
increases  monotonically  with  p.  (The  one  exception  is  that  <t>*(p) 
can  remain  zero  for  two  successive  values  of  p,  in  which  case  the 
second  optimal  point,  x*(p+),  must  be  feasible).  This  is  proven  in 
Lemma  3.3. 

Lemma  3.3  Let  w  :  R  •>  R  satisfy  the  assumptions  of  Lemma  2.2,  and 
let  c^(x),  <t>(x,p),  x*(p),  $*(p)  be  defined  by  Definition  3.1. 

Let  p+  >  p  >  0.  Then 

4>*(P+)  *  **(P) 

with  equality  only  if  c1-(x*(p+))  =  0,  i  =  1,  ....  m. 

Proof:  From  Lemma  2.3 

w(p+-ci(x*(p+)))  >  w(p-ci(x*(p+)))  (3^) 

p-  p 

i  =  1,  ....  m,  with  equality  if  and  only  if  (x*(p+))  =  0  for  all  i. 
Summing  (3.1)  for  i  going  from  1  to  m  yields 
$*(P+)  2  *(x*(p+),p) 


(3.2) 
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with  equality  if  and  only  if  ci(x*(p+))  =  0,  i  =  1,  ....  m. 

Also,  by  the  definition  of  4>*(p) 

#(x*(p+),p)  *  4»*(p)  (3.3) 

Combining  (3.2)  and  (3.3)  completes  the  proof. 

Lemmas  3.2  and  3.3  give  some  insight  into  why  Algorithm  2.3  will 
eventually  solve  the  feasibility  problem  (1.2II).  Since  the  minimum 
value  of  4>(x,p)  is  a  monotonically  increasing  function  of  p,  and 
since  if  4>*(p)  is  positive  for  any  p  then  the  constraints  are  infeasible, 
we  see  that  a  kind  of  "squeeze  play"  is  going  on:  as  p  increases,  either 
<f>*(p)  will  become  positive,  or  if  the  constraints  are  feasible,  <j>*(p) 
will  have  to  stay  nonpositive.  In  the  latter  case,  the  weighting 
function  will  force  x*(p)  to  eventually  become  feasible,  because  it 
will  cause  any  infeasible  point  7  to  eventually  have  <+. (3T,p)  >  0. 

This  argument  is  the  essense  of  our  proofs  that  Algorithm  2.3 
solves  the  feasibility  problem  ( 1 . 2 1 1 ) .  The  remaining  important  fact 
is  that  if  (c.)  is  strictly  feasible  or  strictly  infeasible.  Algorithm 
2.3  will  terminate  satisfactorily  as  soon  as  p  reaches  or  exeeds  some 
finite  value.  This  is  proven  in  Theorems  3.4  and  3.5.  In  Theorem  3.6 
we  prove  a  slightly  weaker  result  in  the  case  when  (c^  is  exactly 
feasible  (feasible  but  not  strictly  feasible).  Recall  that  these 
three  cases  cover  all  possibilities,  as  the  case  of  { c - }  being  infeas¬ 
ible  but  not  strictly  infeasible  is  impossible  as  long  as  D  is  closed 
and  bounded. 

Theorem  3.4  proves  the  finite  termination  of  Algorithm  2.3  if 
{c. }  is  strictly  infeasible.  The  proof  basically  proceeds  by  showing 
that  due  to  the  nature  of  the  weighting  function,  eventually  4. (x ,p)  >  0 
for  any  strictly  infeasible  point  x. 
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Theorem  3.4.  Let  w:  R  -*■  R  be  a  monotonically  increasing  function 
of  y,  and  assume  lim  w(y)  =  l  for  some  finite  l  <  0,  lim  w(y)  =  +  ». 

y-v-co  y-v+oo 

Let  c.(x),  4>(x,p),  x*(p),  <}>*( p)  be  defined  by  Definition  3.1. 

If  (c.j }  is  strictly  infeasible,  then  there  exists  p^  s  0  such  that 
for  any  p  >  p^  4»*(p)  >  0. 

Proof:  Since  { c^ >  is  strictly  infeasible,  there  exists  6^  >  0  such 
that  for  all  x  e  D,  c^ (x)  >  for  at  least  one  i.  Let y^  be  chosen  such  that 
w(y1)  s  -(m-1 )l  (3.4) 

and  define  p^  =  y1  /  Sj.  We  show  that  for  any  p  >  p^,  <}>*( P )  >  0. 

Since  { ci }  is  not  within  6^  of  feasibility,  by  the  monotonicity 
of  w  we  have 

4>*(p)  >  (m-l)f  +  w(p6. ) 

- p - L  <3'5) 

From  the  definition  of  p^  and  y1  and  the  nonotonicity  of  w, 

w ( p 6 x )  >  wCp^p  =  w(yx)  >  -(m-l)f  (3.6) 

Combining  (3.5)  and  (3.6)  completes  the  proof. 

Theorem  3.5  proves  the  finite  termination  of  Algorithm  2.3  if 
{c.j }  is  strictly  feasible.  The  proof  basically  procedes  by  showing 
that  due  to  the  nature  of  the  weighting  function,  eventually 
<j>(x’,p)  >  $(x,p)  for  any  infeasible  point  x  and  any  strictly  feasible 
point  x. 

Theorem  3.5.  Let  w  :  R  ->•  R  be  a  monotonically  increasing  function 
of  y  with  w(0)  =  0,  and  assume  lim  w(y)  =  l  for  some  finite  l  <  0. 

y-*-co 

Let  c.. (x) ,  $(x,p),  x*(p),  $*(,:>)  be  defined  by  Definition  3.1.  If 
{ c^ }  is  strictly  feasible,  then  there  exists  p2  >  0  such  that  for  any 
p  >  p2,  x*(p)  is  feasible. 
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Proof:  Since  (c^ }  is  strictly  feasible,  there  exists  e2  <  0,  x2  e  Rn 

such  that  c-(x2)  s  e2>  1  B  1,  ....  m.  Let  y2  be  chosen  such  that 

w(y2)  <  (m-l)£  (3.7) 

L  m 

and  define  p2  =  y2/e2-  We  show  that  for  a:iy  p  "  p2’  x*^  is  feasible- 
Suppose  that  for  some  p  >  p2,  x*(p)  is  infeasible.  Then  since 
at  least  one  constraint  is  infeasible  at  x*(p)  and  w(y)  >  0  for  all 
y  >  0, 

4>*(p)  >  (m-l)l  .  (3.8) 

P 

Also,  from  the  definition  of  x2  and  the  monotonicity  of  w, 

4>(x2,p)  s  m(w(pe2))  (3-9) 

P 

From  the  definitions  of  p2  and  y2  and  the  monotonicity  of  w, 

w(pej  s  w(p?e?)  =  w(y?)  <  (m-l)l  (3.10) 

c  m 

so  that  from  (3.9)  and  (3.10) 

<t>(x?.p)  <  (m-l)f  (3.11) 

P 

From  (3.8)  and  (3.11),  <f>(x2>p)  <  4>*(p)  which  is  a  contradiction  and 

completes  the  proof. 

The  remaining  case  is  when  { c i :  is  exactly  feasible.  While  this 
case  may  not  be  very  meaningful  in  finite  precision  arithmetic,  results 
indicating  how  Algorithm  2.3  performs  on  it  in  theory  should  indicate 
how  it  will  perform  on  "close  to  exactly  feasible"  problems  in  practice. 
In  Theorem  3.6  we  show  that  if  { c^ }  is  exactly  feasible,  then  given 
any  <5  >  0,  Algorithm  2.3  is  guaranteed  to  produce  an  x*(p)  that  is 
within  6  of  feasibility  as  soon  as  p  reaches  or  exceeds  some  finite 
limit.  We  comment  on  the  optimality  of  this  result  after  the  proof. 

The  proof  of  Theorem  3.6  is  closely  related  to  that  of  Theorem  3.4. 
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Theorem  3.6.  Let  w  :  R  -*■  R  satisfy  the  assumption  of  Theorem  3.5, 

and  assume  in  addition  that  lim  w(y)  =  +  ».  Let  c.(x),  4. ( x »p) 

y  -*•  + 00 

x*(p),  4>*(p)  be  defined  by  Definition  3.1.  If  {ci }  is  feasible, 
then  given  any  <s3  >  0,  there  exists  p^O  such  that  for  any  P^P3» 
x*(p)  is  within  6 3  of  feasibility. 

Proof:  Let  y^  be  defined  by  (3.4)  and  define  P3  =  y We  show 
that  for  any  p  >  p3,  x^(p)  is  within  63  of  feasibility. 

Suppose  that  for  some  p  >  p3,  x*(p)  is  not  within  <$3  of 
feasibility.  Then  c^(x*(p))  >  63  for  at  least  one  i,  and  so  using 
the  monotonicity  of  w, 

♦*(P)  >  (n-l)l  +  w(p  •  63)  (3  12j 

p 

From  the  definitions  of  p3  and  y1  and  the  monotcnicity  of  w, 

w(p  •  63)  >  w(p3  <s3)  =  w(yx)  >  -  (m-l)£  (3.13) 

so  that  from  (3.12)  and  (3.13) 

4>*(P)  >  0.  (3.14) 

On  the  other  hand,  we  have  from  Lemma  3.2  that 
**( P)  *  0 

which  contradicts  (3.14)  and  completes  the  proof. 

Theorem  3.6  implies  that  the  limit  of  the  points  generated  by 
Algorithm  2.3, 

x*  =  lim  x*(p) 

p  ->■  00 

will  be  feasible  for  exactly  feasible  problems.  This  is  a  weaker 
result  than  for  strictly  feasible  problems.  Thus  a  reasonable 
question  to  ask  is  whether  our  result  for  strictly  feasible  problems 
can  be  extended  to  exactly  feasible  problems,  i.e.,  whether  there 
will  exist  any  finite  p  such  that  x*(p)  is  feasible  in  this  case. 
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Example  3.7  shows  that  the  answer  is  no  for  some  exactly  feasible 
problems,  but  yes  for  some  others.  Thus  Theorem  3.6  is  the  best 
result  we  can  prove  about  our  algorithm's  performance  on  exactly 
feasible  problems  in  general,  but  the  algorithm  may  do  better  on 
some  exactly  feasible  problems. 


Example  3.7.  Assume  w  satisfies  Definition  2.1,  let  c^x),  <j>(x,p), 
x*(p)i  4>*(p)  be  defined  by  Definition  3.1,  and  assume  the  set  of 
constraints  {c^  is  exactly  feasible.  Then  it  is  possible  that 
x*(p)  is  infeasible  for  all  p  >  0;  it  is  also  possible  that  x*(p) 
is  feasible  for  all  p  >  0. 

This  is  shown  via  2  examples.  Consider  first 
Cj(x)  =  1  -  x3 
c2(x)  =  X2  -  1 

which  has  its  sole  feasible  point  at  x  =  1.  For  p  =  0,  x*(p)  is 
undefined.  For  p  >  0,  x*(p)  must  satisfy 

p  O  O 

_d_  <j>(x,p)  =  0  =  -3x  w'(p-px  )  +  2xw'(px  -p). 
dx 

It  is  clear  that  x  =  1  is  never  a  solution  to  this  equation,  as  thi 
would  imply  -w^O)  =  0  which  contradicts  the  monotonicity  of  w. 
Thus  x*(p)  must  be  infeasible  for  all  p  >  0.  Incidentally,  for 
sufficiently  large  p. 


x*(p)  «  1 


[ 


w'l°) 


13  p  w  (0)  -  3  w  (0) 


so  that  lim  x*(p)  =  1. 

P  -►  00 


Now  consider 

CjCx)  =  1  -  x3 
c2(x)  =  x3  - 1 
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which  also  has  its  sole  feasible  point  at  x  =  1.  For  p  =  0,  x*(p) 

is  again  undefined.  For  p  >  0 

<j>(x ,p)  =  w(p  •  Cj(x))  +  w(p  •  c2(x)) 

p  (3.15) 

=  w(p  •  ( x ) )  +  w(-p  •  c^ (x) )  . 

p 

From  w(0)  =  0  and  w  strictly  convex,  we  have 
w(a)  +  w( -  a)  _  0 

for  any  a  e  R,  with  equality  if  and  only  if  a  =  0.  Thus  for  any 
p  >  0,  the  unique  minimum  to  (3.15)  occurs  when  p  •  Cj(x)  =  0,  i.e., 
when  c^x)  =  c2(x)  =  0.  Thus  x*(p)  =  1  for  all  p  >  0. 

Theorems  3.4  -  3.6  show  that  Algorithm  2.3  will  solve  the 
feasibility  problem  ( 1 . 21 1 )  correctly.  However,  in  the  exactly 
feasible  case  this  may  only  be  true  in  the  limit  as  p  ■+«>;  correspond¬ 
ingly,  if  { c^ }  is  strictly  feasible  only  for  very  small  (in  absolute 
value)  negative  e,  or  strictly  infeasible  only  for  very  small  positive 
6,  the  correct  solution  may  require  p  to  be  very  large.  In  practice, 
one  probably  wants  to  limit  how  large  p  can  get  in  Algorithm  2.3. 
Therefore,  in  our  computer  implementation  of  this  algorithm,  we  have 
changed  the  termination  criteria  to: 

if  i>*(p)  >  0,  stop  (constraints  are  infeasible) 
if  x*(p)  is  feasible  or  within  5  of  feasibility,  stop  (3.16) 
(constraints  are  considered  feasible) 

where  6  is  a  small  positive  number.  It  is  shown  in  Theorem  3.8  that 
Algorithm  2.3  with  stopping  criteria  (3.16)  is  guaranteed  to  terminate 
as  soon  as  p  reaches  or  exceeds  p  =yj(m+l)/6,  y^  given  by  equation 
(3.4).  The  algorithm  will  be  guaranteed  to  produce  the  correct  result 
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to  feasibility  problem  (1.2II)  unless  {c..}  is  infeasible  but  within 
6  of  feasibility.  If  6  is  a  small  number,  then  for  most  applications 
this  uncertainty  will  be  acceptable,  and  in  some  applications  we  have 
seen  this  "benefit  of  the  doubt"  is  even  desirable.  A  possible  dis¬ 
advantage  is  that  some  applications  really  need  a  feasible  point  if 
one  exists,  and  an  x*(p)  which  is  infeasible  but  within  6  of  feasibility 
will  not  do.  In  the  unlikely  event  that  this  occurs,  one  can  continue 
to  iterate  Algorithm  2.3  in  the  hope  that  it  will  find  a  feasible  point. 

Theorem  3.8.  Let  w  :  R-+R  satisfy  the  assumptions  of  Theorem  3.6, 

and  let  c^x),  $(x,p),  x*(p),  <j>*(p)  be  defined  by  Definition  3.1. 

Then  given  any  <54>0,  there  exists  p4  >0  such  that  for  any  psp4, 

either  <j>*(p)  >  0  or  x*(p)  ‘s  within  64  of  feasibility. 

Proof:  Let  y ^  be  defined  by  (3.4)  and  define  y  =  «4/  (m+1), 

P4  =  yj/y.  Our  proof  is  divided  into  two  cases.  We  show  that  if 

{ c^ }  is  not  within  y  of  feasibility,  then  4>*(p)  >  0  for  any  p  >  p4,  and 

that  if  {Cj}  is  either  feasible  or  within  y  of  feasibility,  that  x*(p) 

is  within  «4  of  feasibility  for  any  p  >  p4- 

If  { }  is  not  within  y  of  feasibility,  then  it  is  shown  in  Theorem3. 4 

that  if>*(p)  >  0  for  any  p  ;>  p4>  Now  consider  the  other  possibility, 

that  { c^ }  is  feasible  or  within  y  of  feasibility.  Then  there  exists 

x4  such  that  c.(x4)  <  y,  i  =  1,  ...,  m,  so  that  for  any  p  >  0, 

4>(x4,p)  s  m  •  w(py)  .  (3.17) 

P 

Now  suppose  that  for  some  p  2  p4,  x*(p)  is  not  within  64  of  feasibility. 
Then  by  the  mono torn’ city  of  w, 

4>*(p)  >  (m-l)£  +  w(p64) 

p 


(3.18) 
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Si  nee  w  is  strictly  convex  and  w(0)  =  0, 

w(p$4)  >  (m+l)w/p<54|  =  (m+1)  w(py)  (3.19) 

m  +1/ 

so  that  combining  (3.18)  and  (3.19) 

$*(p)  >  [(m-l)l  +  w(py)]  +  m « w(py)  .  (3.20) 

P 

Also  by  the  monotonicity  of  w  and  the  definitions  of  and  y^» 

w(py)  s  w(p4y)  =  w(y^)  2  -  (m+l)£  (3.21) 

so  that  combining  (3.20)  and  (3.21), 

4>*( p)  >  m  «  w(Py)  (3.22) 

P 

However,  from  (3.17)  and  (3.22),  <t>(x4»P)<  <f>*(p)  which  is  a  contradic¬ 
tion  and  completes  the  proof. 

A  logical  question  to  ask  now  is:  since  Algorithm  2.3  using 
stopping  criteria  (3.16)  is  guaranteed  to  terminate  for  p^  =  yj(m+l)/6, 
where  y^  is  readily  calculated  given  w  and  m,  why  not  just  solve  the 
linearly  constrained  optimization  problem  (2.4)  with  p  =  p,  instead 
of  executing  Algorithm  2.3  which  solves  (2.4)  for  a  sequence  of  p's? 
There  are  two  reasons.  First,  since  $(x,p)  becomes  more  poorly  behaved 
as  p  increases,  it  may  be  easier  to  solve  (2.4)  for  a  sequence  of  p's 
terminating  with  p  =  p,  using  the  solution  to  each  as  the  starting 
guess  to  the  next,  than  to  just  solve  the  problem  with  p  =  p".  In 
practice,  linearly  constrained  optimization  problems  after  the  first 
in  Algorithm  2.3  require  few  iterations.  Second,  it  is  very  possible 
that  Algorithm  2.3  will  terminate  for  some  p  <  p. 

Another  important  question  which  we  have  not  yet  addressed  is 
whether  our  techniques  can  be  used  to  obtain  a  strictly  feasible 
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point  In  cases  when  {c..}  is  strictly  feasible.  The  answer  is  yes, 
but  how  strictly  feasible  turns  out  to  depend  on  the  weighting  func¬ 
tion  w(y).  This  is  one  aspect  of  the  choice  of  the  weighting  function 
which  is  discussed  in  Section  4. 
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4.  Choice  of  the  weighting  function 

In  Algorithm  2.3  we  have  proposed  a  class  of  algorithms  for 
solving  the  feasibility  problem  (1.2II).  This  class  involves  a 
weighting  function  w(y)  which  must  obey  Definition  2.1,  with  a  specific 
algorithm  involving  a  specific  choice  of  the  weighting  function.  In 
Section  3  we  have  proven  that  any  algorithm  from  this  class  will 
successfully  solve  the  feasibility  problem.  In  this  section,  we  dis¬ 
cuss  the  choice  of  the  weighting  function  w(y).  First,  we  show  that  a 
variety  of  functions  exist  satisfying  Definition  2.1.  Then  we  examine 
which  of  these  is  likely  to  cause  quicker  termination  of  Algorithm  2.3. 
Finally,  we  analyze  (in  Theorem  4.1,  Corollary  4.2  and  Example  4.3) 
the  influence  of  the  weighting  function  on  the  ability  of  Algorithm 
2.3  to  find  an  interior  point  to  a  set  of  strictly  feasible  constraints. 
From  these  considerations,  a  preference  emerges  for  the  weighting  func¬ 
tion  to  use  in  practice. 

The  simplest  (to  write)  function  which  satisfies  Definition  2.1  is 
w^y)  =  ey-  1. 

Indeed,  Wj(y)  and  its  related  multiples  and  powers  are  the  only 

infinitely  differentiable  functions  we  know  of  which  satisfy  Definition 

2.1.  However,  many  twice  continuously  differentiable  functions 

satisfying  Definition  2.1  can  be  constructed  by  "splicing  together" 

two  functions.  For  example,  one  can  let  w(y)  =  y  for  y  <a  <  1,  > 

l-y 

and  then  join  it  to  a  quadratic  for  y  >  o  in  a  way  which  makes  w(y) 
twice  continuously  differentiable.  This  is  done  below  and  drawn  in 
Figure  4.1. 
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w  2(y) 


lL-«f  +  (y-a)  +  _°l_ 

(1-a)3  (1-a)2  1_a 


y  S  a  <  1 


y  >  a 


In  the  same  manner  one  can  form  a  weighting  function  by  splicing 
together  (ey-l)  and  a  quadratic  at  a  point  &  >  0.  This  is  done  in 
w3(y)  and  drawn  in  Figure  4.2.  Function  w3(y)  may  be  preferable  to 
Wj(y)  in  practice  because  it  grows  less  quickly  for  large  y  and  thus 
diminishes  the  chances  of  overflow  or  badly  behaved  optimization 
problems  in  Algorithm  2.3. 
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-  1 


y  <  s  >  0 


w3(y) 


I  6 

lv 


(y-e)‘ 


(y  -  e)  + 


- 1 


y  ^  6 


Figure  4.2 

Certainly  other  functions  satisfying  Definition  2.1  are 
possible,  but  for  our  purposes  it  will  suffice  to  consider  these  three. 
In  particular,  much  of  our  analysis  turns  out  to  concern  only  the  por¬ 
tion  of  w(y)  for  y  <  0,  and  here  the  functions  w.^(y)  =  (/  -  1  and 
w2(y)  =  y/(l-y)  seem  to  be  the  two  important  possibilities. 

One  way  to  compare  the  weighting  functions  is  to  examine  the 
upper  bounds  produced  by  Theorems  3.4  -  3.6  and  3.8  on  how  large  p 
must  be  to  guarantee  termination  of  Algorithm  2.3.  These  bounds  are 
all  functions  of  the  weighting  function,  and  we  are  interested  in 
which  weighting  function  leads  to  the  lowest  bounds.  The  easiest  case 
to  examine  is  the  upper  bound  on  p  when  {c..}  is  strictly  e  feasible 
for  some  e  <  0.  Theorem  3.5  tells  us  that  Algorithm  2.3  is  guaranteed 
to  terminate  for  p  >  y^/e,  where  y^  satisfies 
w(y?)  =  (m  -  1)1  . 

m 
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If  w(y)  =  ey  -  1  for  y  <  0,  it  is  easily  calculated  that  l  =  -  1  and 
y2  =  -  In  m. 

If  w(y)  =  y/(l-y)  for  y  <  0,  again  l  =  -  1  and 
y2  =  -  (m  -1) . 

Since  £nm<m-l  for  all  m  2  2  (and  recalling  e  <  0)  we  see  that  the 
upper  bound  using  e^-l  is  less  than  the  upper  bound  if  one  uses 
y/(l-y).  Note  that  this  analysis  only  involves  the  portion  of  w(y) 
for  y  <  0  .  Thus  using  w(y)  =  ey  -  1  for  y  <  0  may  be  preferable  to 
using  w(y)  =  y/(l-y)  for  y  <  0  on  strictly  feasible  problems. 

Now  consider  the  bound  given  by  Theorem  3.4  in  the  case  when 
{c.j}  Is  strictly  6  infeasible  for  some  6  >  0.  The  theorem  says  that 
Algorithm  2.3  is  guaranteed  to  terminate  for  p  2  y^/  &,  where  y ^ 
satisfies 

w(yx)  =  -(m-l)£. 

If  w(y)  =  ey  -  1  then 

y1  =  £nm. 

If  w(y)  =  w2(y)  then 


r 


m  -  1 


yl 


3  a  -  1  +  1  -  a  V4  m  ( 1-ci)  -  3 

2  2 


a  2  m  -  1 
m 


a  s  m  -  1 
m 


Thus  the  situation  here  is  less  clear  cut  than  the  above.  If 
a  >  (m-l)/m,  then  the  upper  bound  for  w2(y)  clearly  is  less  than 
the  upper  bound  for  Wj(y).  If  a  <  (m-l)/m  then  the  bound  for  w2(y) 
may  exceed  the  bound  for  w^(y).  The  bound  for  w^(y)  is  always  greater 
than  or  equal  to  the  bound  for  w^(y).  Note  that  this  analysis  only 


involves  w(y)  for  y  >  0. 
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However,  the  selection  of  the  weighting  function  for  y  >  0  is 
complicated  by  the  presense  of  two  conflicting  considerations:  the 
more  quickly  w(y)  increases,  the  lower  will  be  the  upper  bound  on  p 
required  by  Theorem  3.4,  but  the  larger  will  be  the  chance  for  over¬ 
flow  or  bad  behavior  in  the  objective  function  $(x,p)  if  p  does  be¬ 
come  large.  So  while  we  could  always  set  a  >  (m-1)  /m  and  thus  have 
a  smaller  y^  for  w^(y)  than  for  w^(y)  or  w^y),  we  might  not  want  to 
do  this  because  it  might  make  w(y)  too  steep  for  y  >  0.  Similarly, 
Theorems  3.6  and  3.8  show  that  if  we  stop  our  algorithm  when  x*(p) 
is  within  6  of  feasibility,  then  the  upper  bound  on  p  is  again  pro¬ 
portional  to  y  ,  but  the  same  conflicting  considerations  cloud  the 
selection  of  w(y).  Thus  we  believe  that  only  experimentation  can 
determine  which  weighting  function  is  preferable  for  y  >  0. 

The  other  important  regard  in  which  the  selection  of  the  weight¬ 
ing  function  can  influence  the  performance  of  Algorithm  2.3  is  the 
ability  of  the  algorithm  to  find  strictly  feasible  points  when  { c^. } 
is  strictly  feasible.  While  Algorithm  2.3  as  written  will  stop  as 
soon  as  a  feasible  point  is  found,  some  applications  would  like  a 
point  which  is  as  feasible  as  possible.  Therefore,  we  might  be 
curious  whether,  if  we  continue  to  run  Algorithm  2.3  after  it  has 
found  a  feasible  point  to  a  strictly  t  feasible  problem,  it  will 
eventually  find  an  x*(p)  which  is  "almost  as  feasible  as  possible," 
i.e.,  strictly  i  feasible  for  *  greater  than  but  arbitrarily  close 
to  e .  Theorem  4.1  enables  us  to  answer  this  question,  by  giving  a 
necessary  condition  for  Algorithm  2.3  to  eventually  produce  a  strictly 
c  feasible  x*(p)  on  a  strictly  e  feasible  problem  for  a  given  c  in  (r  ,o). 
This  condition  involves  w(y)  only  for  y  <  0.  Corollary  4.2  then  shows 
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that  w(y)  =  ey  -  1  satisfies  this  condition  for  any  e  in  (e,0).  Thus 
Algorithm  2.3  with  w(y)  =  e^  -  1  for  y  <  0  will  produce  points  which 
are  "almost  as  feasible  as  possible."  However,  w(y)  =  y/  (1-y)  does 
not  satisfy  the  condition  of  Theorem  4.1  for  e  <  e/m,  and  indeed  it 
is  shown  in  Example  4.3  that  continuing  Algorithm  4.3  with 
w (y)  =  y/  (1-y)  for  y  <  0  may  fail  to  produce  points  which  are 
"almost  as  feasible  as  possible." 

Theorem  4.1 ,  Let  w  obey  the. assumptions  of  Theorem  3.5  and  define 

w(y)  =  w(y)-£.  Let  c^x),  <p(x,p),  x*(p),  <>*(p)  be  defined  by 

Definition  3.1.  Let  e  <  e  <  0,  and  assume  (c^x)}  is  strictly  e 

feasible.  Suppose  there  exists  p,.  >  0  such  that  for  all  p  >  p^, 

m  •  w(pe)  <  w(pc) .  (4.1) 

Then  for  any  p  >  p5,  x*(p)  is  strictly  e  feasible. 

Proof :  Suppose  there  exists  p5  satisfying  the  conditions  of  the 

theorem,  and  that  for  some  p  >  p^,  x*(p)  is  not  e  feasible.  Then  at 

least  one  constraint  has  value  greater  than  e  at  x*(p),  so  that 

$*(p)  >  (m-1)  l  +  w(p  e)  =  m  •  l  +  w(p  e  )  (4.2) 

P  P 

However,  since  {c..}  is  strictly  e  feasible,  there  exists  x^  e  Rn 

such  that  c.(xc)  <  e,  i  =  1,  ...,  m.  Thus  from  the  definition  of 

x*(p)  and  the  monotonicity  of  w, 

4>*(p)  s  mw(p  e)  =  m  »  l  +  mw(p  e)  (4.3) 

P  P 

Equations  (4.2)  and  (4.3)  imply  that 
w(p  t)  <  m  •  w(p  c) 

which  contradicts  (4.1)  and  completes  the  proof. 

Corollary  4.2.  Let  w(y)  =  ey  -  1  for  y  <  0,  and  let  c^x),  <j>(x,p). 
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x*(p)  be  defined  by  Definition  3.1.  Let  e  <  c  <  0  and  assume  that 
{ }  is  strictly  e  feasible.  Then  there  exists  Pg  s  0  such  that  for 
all  p  >  Pg,  x*(p)  is  strictly  e  feasible. 

Proof;  From  Theorem  4.1,  it  suffices  to  show  that  there  exists  pg 
such  that  for  all  p  >  Pg,  (4.1)  is  true.  For  w(y)  =  e^-l,  w(y)  =  ey, 
so  (4.1)  requires  that 

me  p  e  <  ep  " 


or 


sP(e  '  e) 


m  s  er 

This  clearly  is  true  for  any  p  >  In m  =  pK. 

$  o 


k  -  e 


Now  let  us  consider  what  Theorem  4.1  says  about  the  weighting 
function  w(y)  =  y/(l-y),  y<0.  In  this  case  w(y)  =  l/(l-y),  and  so 
(4.1)  requires 

m  <  1 

1-pe  Re 


or 


m-1  s  p(m  e  -  e) 

which  is  only  possible  if 

t  >  E_  . 

m 

Thus  with  this  weighting  function.  Theorem  4.1  only  guarantees  that 
continuing  to  iterate  Algorithm  2.3  on  a  strictly  e  feasible  problem 
will  eventually  produce  an  x*(p)  which  is  strictly  e  feasible  for 
e  >  e/m.  Indeed,  in  Example  4.3  we  give  an  example  where  Algorithm 
2.3  does  not  produce  any  x*(p)  which  is  strictly  (-2.5)  feasible  on 
a  problem  which  is  strictly  (-3)  feasible.  For  simplicity  we  let  the 
constraints  be  linear;  more  extreme  examples  are  possible  for 
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7  +  24p _ 1 

J~ 7  +  24p)2  -  (6  +  18p)2  J 

which  is  clearly  greater  than  zero  for  all  p  >  0.  finally 

litn  h(p)  =2-  4/7  =  .488  =  x*  . 

p-»-“  7 

Since  c^(x*)  =  -2.488,  it  is  clear  from  Figure  4  that  no  x^(p) 
is  strictly  e  feasible  for  any  e  <  -2.5. 

Thus  with  regard  to  finding  interior  points,  the  choice  w(y)  = 
ey  -  1  for  y  <  0  seems  superior  to  w(y)  =  y/  (1-y).  This  concurs  with 
our  analysis  earlier  in  this  section.  For  y  >  0  the  issue  "'s  unsettled 
by  theoretical  analysis.  Thus  in  practice  one  probably  wants  to 
experiment  between  w^Cy),  w„(y),  and  possibly  a  hybrid  which  com¬ 
bines  ey-  1  for  y  <  0  with  a  function  like  w2(y)  for  y  >  0. 

It  should  be  mentioned  that  the  class  of  weighting  functions 
satisfying  Definition  2.1  is  not  the  broadest  class  of  twice  con¬ 
tinuously  differentiable  functions  for  which  the  results  of  this 
paper  are  true.  For  example,  I.  Zang  [13]  has  suggested  to  us  that 
we  consider  weighting  functions  which  obey  Definition  2.1  except 
that  they  are  linear  in  some  intervals.  An  example  (related  to  one 
in  [14])  is  given  below  and  drawn  in  Figure  4.4. 


h'(p)  = 


3p‘ 


Figure  4.4 


The  class  of  functions  which  look  like  this,  that  is,  obey  Definition 
2.1  except  that  they  may  be  strictly  convex  only  in  some  open  neighbor¬ 
hood  of  y  =  0  and  convex  elsewhere,  is  a  class  to  which  all  the  results 
of  this  paper  are  readily  extendable.  It  may  prove  to  be  attractive 
computationally. 
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5.  Computational  results 

In  this  section  we  discuss  some  important  aspects  of  the 
computer  implementation  of  our  algorithm  for  solving  feasibility 
problem  (1.2II),  and  present  some  sample  computational  results.  Due 
to  the  length  of  this  paper,  and  because  we  are  continuing  to  revise 
and  experiment  with  the  software  for  our  algorithm,  our  treatment  in 
this  section  is  brief.  A  future  paper  will  be  specifically  devoted 
to  the  implementation  and  test  results  of  our  algorithm. 

In  implementing  Algorithm  2.3,  three  important  aspects  need  to 
be  considered: 

1.  the  choice  of  the  weighting  function 

2.  the  linearly  constrained  optimization  algorithm,  especially 
the  strategy  to  attempt  to  find  the  global  minimum 

3.  the  strategy  for  increasing  p. 

The  choice  of  the  weighting  function  was  discussed  in  Section  4.  In 
our  implementation  so  far  we  have  used  w(y)  =  ey-l,  which  is  supported 
by  our  analysis  there.  We  have  not  had  any  trouble  with  overflows 
or  badly  conditioned  optimization  problems.  As  is  mentioned  in 
Section  4,  we  also  intend  to  experiment  with  the  exponential -quadratic 
hybrid  w-j(y). 

The  linearly  constrained  optimization  algorithm  is  the  heart  of 
our  computer  algorithm.  Our  requirements  for  this  algorithm  are  dif¬ 
ferent  enough  from  the  standard  ones  that  we  deemed  it  preferable  to 
write  our  own  code  rather  than  modify  an  existing  one.  In  the  initial 
version  of  our  algorithm,  we  have  developed  an  algorithm  similar  to 
the  ones  proposed  by  Fletcher  [6]  and  Goldfarb  [8].  The  algorithm 
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uses  finite  difference  approximations  to  the  first  and  second  deriva¬ 
tives  of  $(x,p)  (i.e.,  of  c^(x),  i  =  1,  ....  m)  because  this  is  reason¬ 
able  in  the  applications  with  which  we  are  familiar.  If  the  projected 
Newton  direction  is  not  a  descent  direction  then  it  takes  a  step  in 
the  projected  steepest  descent  direction,  although  clearly  this  aspect 
of  the  algorithm  can  be  improved.  Otherwise,  it  attempts  to  find  a 
local  minimum  of  $>(x,p)  with  respect  to  x,  in  a  standard  manner.  Along 
the  way,  it  attempts  to  obtain  information  indicating  whether  a  lower 
value  of  (|>(x,p)  might  exist  elsewhere.  Information  is  gathered  from 
directions  of  negative  curvature  (if  any  are  encountered)  indicating 
regions  of  the  domain  in  which  a  lower  value  might  be  found.  Then, 
after  a  local  minimum  has  been  found,  this  information  together  with 
a  sampling  of  points  near  and  far  from  the  local  minimum  is  used  to 
determine  whether  a  restart  of  the  algorithm  in  another  region  is 
likely  to  produce  a  lower  value  of  <>(x,p).  If  so,  the  algorithm  is 
restarted  and  this  entire  process  is  repeated  (possibly  more  than 
once).  Such  a  procedure  is  certainly  ad-hoc  and  will  not  work  on 
every  problem,  but  it  has  proven  successful  on  the  majority  of  problems 
with  multiple  local  minima  that  we  have  tried. 

Our  linearly  constrained  optimization  algorithm  also  requires 
some  special  features  due  to  the  setting  in  which  it  is  used.  Certainly, 
if  at  any  stage  in  the  minimization  of  <j>(x,p)  the  current  iterate  x£ 
is  feasible,  then  there  is  no  need  to  continue  the  minimization  and 
so  both  the  linearly  constrained  optimization  algorithm  and  Algorithm 
2.3  are  terminated.  Secondly,  while  the  value  of  4>(x,p)  is  bounded 
below  for  any  p  >  0,  it  is  possible  that  it  has  no  finite  minimum 
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point  x*  if  the  domain  D  is  not  bounded.  This  is  treated  in  the  standard 

way  by  terminating  the  local  minimization  of  <{>(x,p)  if  successive 

iterates  of  x  lead  to  insufficient  decrease  in  <j>(x,p).  Finally,  when 

m 

p  =  0  it  is  possible  that  <j>(x ,0)  =  £  c.(x)  is  unbounded  below  if  D 

i=l  1 

is  unbounded,  and  that  the  algorithm  may  pursue  an  unbounded  path  with¬ 
out  considering  any  feasible  points.  Two  simple  examples  are 
Cj(x)  =  x2-  1 
c2(x)  =  x3 

which  is  feasible  and 

c^x)  =  x2  +  1 
c2(x)  =  x3 

which  is  infeasible.  In  both  cases  4>(x .0)  decreases  without  limit 
as  x  -*■-<*>  ,  and  the  algorithm  will  pursue  this  path  without  locating  a 
feasible  point  if  the  initial  x  is  less  than-1.  Thus  our  linearly  con¬ 
strained  optimization  algorithm  treats  the  case  p  =  0  specially,  and 
in  particular,  once  ^(x,0)  <  0,  it  detects  situations  like  the  above 
and  terminates  the  minimization  of  <}>(x,0)  prematurely.  Algorithm  2.3 
then  increases  p  and  continues. 

Our  entire  algorithm  for  linearly  constrained  optimization  in  the 
context  of  Algorithm  2.3  is  currently  undergoing  major  revision,  and 
will  be  reported  in  more  detail  on  in  the  future.  Our  initial  algorithm 
has  worked  well  enough  for  our  implementation  of  Algorithm  2.3  to 
be  quite  successful. 

The  strategy  for  increasing  p  has  been  given  considerable 
attention,  although  perhaps  surprisingly  Algorithm  2.3  does  not  seem 
too  sensitive  to  it.  At  each  iteration  of  Algorithm  2.3,  we  model  the 
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function  <j>*(p)  =  $(x*(p) ,p)  by  a  function  m(p),  and  then  choose  the  next 
value  p+  such  that  m(p+)  =  0.  In  this  way  we  attempt  to  "force  the 
issue"  by  either  producing  a  value  of  $*(p)  >  0,  or  a  feasible  x*(p). 

The  model  function  m(p)  is  designed  to  fit  (5.1) 

m(pc)  =  *+(Pc).  <"(pprev)  =  ^Pprev*’  m'(Pc>  =  <_♦ Wp)  .p)  |p=Pc 

where  pc  and  Pprev  are  the  current  and  immediately  previous  values 
of  p.  The  derivative  in  (5.1)  is  an  easily  obtainable  approximation  to 
4>*'(pc);  in  a  number  of  tests  taking  finite  difference  approximations 
to  the  latter,  we  found  the  two  derivatives  to  be  almost  the  same.  While 
one  can  fit  (5.1)  with  a  quadratic,  it  seems  preferable  to  make  use  of 
the  problem  structure.  We  currently  use 
m(p)  =  aw  (p  •  e)  +  y 

where  a,  6  and  y  are  chosen  to  satisfy  (5.1);  a  perhaps  more  attractive, 

but  more  difficult  to  obtain,  possibility  is 

m(p)  =  aw(p  •  B)  +  y 
P 

If  the  predicted  p+  is  not  in  the  range  [2pc>  10pc]  we  use  aPP'"o- 
priate  extremum  of  this  range  as  p+.  in  the  initial  case  (p£  =  0) 
the  model  is  constructed  somewhat  differently. 

In  Tables  5.1  -  5.4  we  present  four  indicative  runs  of  our 
initial  implementation  of  Algorithm  2.3.  Mainly  they  illustrate  that 
the  algorithm  works  in  practice  about  as  is  predicted  by  our  theory. 

What  is  perhaps  most  striking  in  our  tests  so  far  is  how  few  values 
of  p  seem  to  be  required  to  obtain  the  correct  answer.  Table  5.1 
and  5.2  illustrate  the  algorithm  on  feasible  problems.  Table  5.1  is 
a  case  where  the  global  characteristics  of  our  linearly  constrained 


optimization  were  required  due  to  the  shape  of  the  tilted  sine  function, 
Tables  5.-3  and  5.4  illustrate  the  algorithm  on  infeasible  problems. 
Table  5.4  is  interesting  because  the  optimal  point  x*(p)  remains  the 
same  for  all  p,  but  p  must  reach  a  positive  threshhold  value  before 
4>*(p)  >  0  and  infeasibility  is  detected. 


Table  5.1 


Cj(x)  =  sin(xx2  +  x22)  +  x1 
C£(x)  =  — £  (xj+^r)2  +  x22  -1 

IT 

initial  guess  =  (0,  (3n/2)1//2) 


p 

iterations  inside 

X*(P) 

<Mp) 

linearly  constrained 

optimization  routine 

0 

14 

(-  6.3,  0.33) 

-  5.2 

.37 

2 

(-  6.3,  0.01) 

-  2.3 

.7 

1 

(-  6.0,  0.C1) 

-0.47 

(feasible  point  found) 


constraint  values 
at  x*(p) 


-  5.3,  0.1 
-5.8, 0.1 
-6.8,  -0.36 


Table  5.2 


Cj(x) 

= 

(Xj/2  -  3)  Xj  +  2x2 

-  1 

c2(x) 

=  X1 

+  (x2/2  -  3)  x2  +  2x3 

-  1 

c3(x) 

=  x2 

+  (x3/2  -  3)  x3  +  2x4 

-  1 

c4(x) 

=  X3 

+  (*4/2  -3)  x4  +  2xg 

-  1 

c5(x) 

=  x4 

+  (x5/2-3)  x5 

-  1 

initial  guess 


(1.1. 1.1.1) 
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iterations  inside 

x*(?) 

$*(p) 

constraint  values 

linearly  constrained 

at  x*(p) 

optimization  routine 

1 

(2,  0,0,0,  1) 

-7.5 

-5,1-1 ,1,-3. f; 

2 

(1 .3,0.7 ,0.3, 0.2, 0.3) 

-4.9 

-2. 5, -1,-0. 7, 

-0.6, -1.7 

( feasible  point  found) 

Table  5.3 

2  2 

c^x)  =  Xj  +  2x2  -  4 

?  3 

C2(x)  =  x^  +  2x2  +  x^  -  8 

c3(x)  =  (x^l)2  +  (2x2-21/2)2  +  (x3-5)2  -4 

initial  guess  =  (1.0, 0.7, 5.0) 


P 

iterations  inside 

X*(P) 

Mp) 

constraint  values 

linearly  constrained 

at  x*(p) 

optimization  routine 

0 

4 

(0.3, 0.3, 1.5) 

1.7 

-3.7,  -3.7,  9.1 

(concluded  infeasible) 
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10 


Table  5.4 

c^(x)  =  quadratic  with  minimum  value  of  1  at  x  =  0  e  R 
c2(x),  ....  c1Q(x)  =  quadratic  with  minimum  value  of  -1  at  x  =  0  e  f 
initial  guess  =  (1,1, ...1) 


p 

iterations  inside 

X*(P) 

4>*(p) 

linearly  constrained 

optimization  routine 

0 

1 

*.(0,0, ...0) 

-9.0 

0.4 

1 

o 

o 

o 

<> 

-3.7 

4.0 

1 

o' 

* 

o 

o 

11.2 

constraint  values 
at  x*(p) 

(l.-l.-l . -1) 

(l.-l.-l . -1) 

(1,-1, -1.....-1) 


10 
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6.  Summary  and  directions  for  further  work 

We  have  presented  a  new  class  of  algorithms  (Algorithm  2.3)  for 
finding  a  feasible  point  or  determining  infeasibility  of  a  set  of 
nonlinear  inequality  constraints 

c^x)  <0,  i  =  1,  ....  m 

subject  to  satisfying  a  set  of  linear  inequality  constraints 
£j(x)  s  0,  j  =  1,  ....  k 

which  are  known  to  be  feasible  (problem  (1.2II)).  Our  class  of 
algorithms  consists  of  solving  a  series  of  problems  of  the  form 

1  m 

min  4>(x,p)  «  -  l  w(p*c.(x))  . 

xeRn  p  i-1  1  '0,i 

subject  to  £.(x)  s  0,  j  =  1,  . .. ,  k 

J 

for  increasing  values  of  the  nonnegative  parameter  p,  where  w  is  a 
weighting  function  obeying  Definition  2.1.  We  have  shown  that  in 
almost  all  cases  (the  lone  exception  being  a  problem  which  is  exactly 
feasible,  i.e.,  feasible  but  with  no  point  where  all  the  constraints 
are  negative)  that  once  p  reaches  or  exceeds  a  finite  value,  the  answer 
to  (6.1)  will  either  be  a  feasible  point  or  establish  infeasibility  of 
our  problem.  This  may  depend,  however,  on  the  ability  of  our  algorithm 
for  (6.1)  to  find  the  global  solution  to  this  problem.  An  initial 
computer  implementation  of  our  algorithm  using  w(y)  =  -  1  has 

successfully  solved  a  wide  range  of  feasible  and  infeasible  problems. 

We  are  actively  continuing  work  on  many  aspects  of  this  project. 
Foremost  is  the  improvement  of  the  computer  implementation  of  Algorithm 
2.3.  The  main  part  of  this  job  is  the  continued  development  of  a 
linearly  constrained  optimization  routine  which  works  on  nonconvex  as 
well  as  convex  problems  and  attempts  to  find  the  global  minimum.  We 
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also  intend  to  experiment  with  other  choices  of  the  weighting  function 
w(y)  as  indicated  in  Section  4,  and  with  other  strategies  for  increas¬ 
ing  the  parameter  p  including  the  final  one  indicated  in  Section  13. 
Ultimately  we  hope  to  develop  reliable  and  extensively  tested  soft¬ 
ware  for  solving  problem  (1.2II). 

A  second  remaining  aspect  of  our  work  is  to  compare  our  algorithm 
with  other  approaches  that  can  be  used  to  solve  the  feasibility  problem. 
In  particular,  we  have  in  mind  approaches  involving  non-differentiable 
optimization  (e.g.,  the  minimax  formulation  of  equation  (2.3))  or 
optimization  subject  to  nonlinear  constraints  (e.g.,  the  fo.mulation  of 
Bertsekas[  1]).  We  would  naturally  be  interested  in  comparison  with 
any  other  algorithm  which  solves  the  feasibility  problem. 

Finally,  we  would  like  to  extend  the  work  of  this  paper  to 
cover  a  broader  class  of  feasibility  problems.  In  particular,  we 
have  excluded  equality  constraints  from  our  consideration  so  far. 

While  the  addition  of  linear  equality  constraints  causes  no  change 
to  our  theory  and  very  little  to  our  software,  how  to  handle  nonlinear 
equality  constraints  is  far  less  obvious.  Certainly  a  nonlinear 
equality  constraint  c^(x)  =  0  can  be  equivalently  expressed  as  the 
two  inequality  constraints 

cl(x)s0  (6.2) 

-  c^x)  <  0  . 

However,  a  system  of  inequality  constraints  including  (6.2)  is 
exactly  the  type  which  may  give  our  algorithm  trouble,  because  any 
feasible  point  will  have  to  be  exactly  feasible.  The  situation  may 
not  be  as  bad  as  this;  in  particular,  a  generalization  of  the  second 
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part  of  Example  3,7  shows  that  if  a  system  of  inequality  constraints 
contains  only  pairs  of  the  type  (6.2)  and  is  feasible,  then  the  x*(p) 
produced  by  Algorithm  2.3  will  be  feasible  for  any  p  >  0.  However, 
if  standard  inequality  constraints  are  included  as  well,  then  x*(p) 
may  be  feasible  only  in  the  limit  as  p-»-°°.  Of  course,  if  we  stop 
our  algorithm  when  4>*(p)  >  0  or  x*(p)  is  within  6  of  feasibility 
for  some  6  >  0  (stopping  criteria  (3.16))  then  it  will  terminate  for 
finite  p  on  these  problems  as  well,  and  the  practical  question  is 
whether  p  has  to  grow  too  large.  It  will  be  interesting  to  see  how 
our  computer  algorithm  performs  on  systems  with  equality  constraints 
expressed  in  form  (6.2),  and  also  to  see  whether  any  other  approaches 
lead  to  a  better  treatment  of  equality  constraints. 
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